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1.0  INTRODUCTION 

Analysis  of  gas  flows  in  which  there  is  significant  energy  transfer  by  radiation  requires 
a  method  of  analysis  of  the  radiant  energy  transfer  field.  Analytical  solutions  of  the  integral 
equation  of  radiant  energy  transport  have  so  far  been  limited  to  a  very  few  simple  geometries, 
one  of  which  is  the  field  with  cylindrical  symmetry  —  a  good  match  with  the  geometry  of 
modern  high-pressure  arc  heaters,  provided  the  electric  discharge  remains  on  the  centerline 
of  the  heater.  Unfortunately,  the  arc  column  in  such  devices  undergoes  substantial  excursions 
from  the  centerline  a  large  part  of  the  time.  It  is  therefore  desirable  to  develop  analytical 
methods  capable  of  dealing  with  strong  temperature  asymmetries.  As  part  of  more  general 
flow-field  solutions,  development  of  a  means  of  analysis  of  nonsymmetric  radiation  fields 
is  required. 


2.0  RADIATIVE  ENERGY  TRANSFER  IN  HOT  GASES 


Analysis  of  behavior  of  very  hot  gases  requires  the  modification  of  the  equations  of  motion 
by  the  addition  of  a  radiant  energy  transfer  term  to  the  conservation  of  energy.  Transport 
of  radiant  energy  in  a  region  of  hot  gas  is  described  by  the  spatial  and  directional  variation 
of  the  spectral  radiation  intensity  function,  I,,  (r,  T).  This  function  is  defined  as  the  rate  of 
flow  of  radiant  energy  of  frequency  v  in  the  direction  of  the  unit  vector  f  at  a  point 
specified  by  r,  per  unit  area  normal  to  the  direction  of  propagation,  per  unit  solid  angle 
increment  about  the  direction  of  propagation,  and  per  unit  increment  of  photon  frequency 
about  p.  The  spatial  variation  of  1^  in  a  given  direction  is  described  by  the  radiative  transfer 
equation  (RTE)  for  a  nonscattering  gas: 

-^  =  -  K,  [I,  -  Bj  (1) 

where  B„  is  the  equilibrium  (blackbody)  intensity,  K„  is  an  attenuation  or  absorption 
coefficient,  and  S  is  the  distance  measured  in  the  direction  of  propagation. 


The  “gray  gas”  approximation  is  frequently  justifled;  namely,  the  absorption  coefficient 
is  independent  of  photon  frequency,  K,.  =  K;  hence,  the  intensities  can  be  integrated  over 
all  frequencies  to  become: 


I  «  j"  I,  dp  and  B  s  j*  B,  d.-  =  a  TV*- 


With  this  approximation,  the  RTE  has  the  form: 


dl 

dS 


(2) 


(3) 
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The  radiant  intensity  at  any  point  in  any  given  direction  of  propagation  is  obtained  by 
integrating  this  ordinary  differential  equation  along  that  direction  of  propagation  from  the 
boundary  of  the  radiating  gas  to  the  (observation)  point,  represented  formally  by: 

1^.0  =  I(rB)e-^  +  -  JJ®T4e--dT  (4) 

IT  * 


where  the  variable  r  is  the  (dimensionless)  optical  depth  along  the  inverse  speciHed  direction 
of  propagation,  referenced  to  the  observation  point  (see  Fig.  la) 

r.|®K<b  (5) 

The  two  terms  on  the  right  in  the  formal  solution,  £q.  (4),  have  the  following  interpretation. 
The  first  term  is  the  radiation  intensity  in  direction  of  f  at  a  boundary  point  lying  on  the 
inverse  projection  of  f  through  the  observation  point,  attenuated  by  absorption  as  measured 
by  the  optical  thickness  between  the  two  points.  The  second  term  is  the  summation  of 
contribution  of  emission  in  the  direction  of  f  from  all  points  lying  between  the  boundary 
and  observation  points,  each  elementary  contribution  attenuated  by  absorption  attributable 
to  the  corresponding  optical  thickness  between  each  point  and  the  observation  point. 


The  radiant  energy  (heat)  flux  vector  in  the  radiating  gas  is  given  by  an  integration  of 
the  intensity  function  over  all  possible  directions: 


q|‘=  jj’ I  (f,  £ )  i;  dO 


(6) 


where  ^  is  a  direction  cosine  of  the  unit  vector  £ .  The  radiant  energy  transport  term  required 
for  the  energy  equation  is  the  divergence  of  the  vector  q^: 


V .  ^  j*- 1,  do  =  ( 

^  dxi  axi  )o  ‘  J 


from  Eq.  (3), 


dxi 
V  •  q^ 


dl 


=  -KjiJ’irr. 


“  dS 

I 

n  dQ  -  4aT^ 


dQ  and. 


(7) 


The  radiant  intensity  in  both  Eqs.  (6)  and  (7)  is  obtained  from  the  formal  solution  of  the 
RTE,  Eq.  (4): 


’4t 

0 

I(TB)e-"B  +  -  1?  e-"  dr 

TT 

tiiU 

V  .  qR  =  -K 

ir  j^T^e-'dr 

dO  -  4aT* 

(8a) 

(8b) 
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Examination  of  these  relations  shows  that  the  extreme  complexity  in  calculating  or 
V  •  qK  for  any  but  the  simplest  geometry  lies  in  obtaining  integral  solutions  of  the  RTE  for 
a  very  large  number  of  directions  of  propagation  at  each  point,  to  provide  a  well-deHned 
integrand  for  the  integrals.  (In  principle,  the  RTE  must  be  satisfied  in  all  directions  in 
a  Ak  solid  angle  at  each  point).  This  results  from  the  fact  that  the  radiant  intensity  is  not 
a  vector  quantity  definable  in  terms  of  only  three  components. 

3.0  EXACT  ANALYTICAL  SOLUTIONS 

Clearly,  Eq.  (8b}  forms  a  very  complicated  integro-differential  equation  for  the  radiant 
heat  addition  term  in  the  energy  conservation  equation.  Exact  analytical  solutions  have  been 
limited  to  very  simple  cases  involving  a  single  independent  space  variable,  either  one- 
dimensional  or  axisymmetric  cases.  (“Exact**  solution  here  means  that,  within  the  gray  gas 
approximation,  Eqs.  (Ba)  and  (8b)  are  satisfied  for  all  directions  of  propagation  at  all  points 
of  the  radiation  field). 

3.1  ONE-DlMENSlONAL  (PLANAR)  CASE 

The  radiant  field  in  a  gray  gas  lying  between  two  plane  parallel  black  walls  of  infinite 
extent  has  been  subject  to  exact  analysis  (Refs.  1 — 2).  The  location  between  the  two  plates 
is  defined  by  the  optical  thickness  function  (see  Fig.  lb), 

i,{x)=|jKdx'  (9) 

with  X  measured  normal  to  the  wall.  Note  the  difference  between  this  function  and  the  optical 
depth  function,  Eq.  (S),  which  is  evaluated  in  the  direction  of  propagation.  The  radiant  heat 
addition  term  in  this  case  is  found  to  be: 


where  subscripts  R  and  L  refer  to  right  and  left  walls,  i; '  is  a  running  variable  of  integration 
normal  to  the  wall,  f  is  the  cosine  of  the  angle  between  the  propagation  direction  and  the 
normal  to  the  walls,  and  T  =  T(ii).  Integrations  over  i  represent  the  integration  of  intensity 
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over  all  possible  directions  in  a  4t  solid  angle.  These  integrals  belong  to  the  family  of  integro- 
exponential  functions  defined  by: 


En(Z)  =  [‘e.-2e-z/rdf 


(11) 


which  functions  have  been  evaluated  and  tabulated  by  various  authors  (Ref.  3).  In  terms 
of  these  functions,  the  heat  addition  term  is: 


-  =  2  Ka  E2  in)  +  jj  T<  E|  (ij-if')  dij'  +  Ej  -  i,)  (12) 


+  E,  (,'  -  ,)  d,'  -  2T4] 

The  radiant  heat  flux  vector  itself  is  normal  to  the  walls  with  a  magnitude: 

qR  =  2  <r  [Ti,L  E3  (ij)  +  jj  T*  Ej  (?y  -  ij')  dif'  -  E3  (ij^r  -  ij)  - 

(13) 

|’'*'“T4E2(V  -  ,)di,'] 

As  noted,  Eqs.  (12)  and  (13)  provide  the  exact  solution  of  the  radiative  transfer  equation 
for  the  plane-wall  case,  with  the  gray  gas  approximation.  A  family  of  such  solutions  for 
radiative  heat  flux  between  two  parallel  walls  is  plotted  in  Fig.  2  for  a  range  of  optical  thickness 
from  0.244  to  12.0  and  for  a  linear  gas  temperature  distribution  from  25,000  to  5,000“R. 
This  range  of  optical  thickness  was  obtained  by  holding  the  absorption  coefficient  constant 
(K  =  24.4  ft“'  =  0.80052  cm"*)  and  allowing  wall  spacing  (A),  to  vary  from  0.12  to  6.0 
in.  Inspection  of  Eq.  (13)  shows  that  zero  optical  thickness  would  give  a  constant  heat  flux, 
q*t  =  a  [Twl^  —  Twr^J  =  185,500  Btu/ft^sec  at  all  points  between  the  walls.  As  the  optical 
thickness  is  increased,  the  heat  flux  decreases  because  of  absorption,  and  a  flux  maximum 
is  established  at  some  point  near  the  hotter  wall. 

In  a  frequently  used  approximation,  the  integro-exponendals,  Eq.  (11),  are  replaced  by 
simple  exponential  functions: 

El  (Z)  s  abe-oz  Ej  (Z)  a  ae-w,  E3  (Z)  s  A  e-*>z  (14) 

b 

with  a  =  3/4,  b  =  3/2.  The  basic  goal  of  an  “exact”  solution,  namely,  the  satisfaction 
of  the  RTE  in  all  possible  directions,  is  not  compromised  to  any  substantial  degree  by  use 
of  this  approximation.  The  heat  flux  expression  becomes: 
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2  a 


T^l  y  +  7  l;  e-5/2<,-v)  d„' 

_  A  j’WR  j4  g-3/2l,'-,) 


(15) 


The  radiant  heat  flux  between  two  parallel  plane  walls  obtained  from  Eq.  (15)  is  compared 
with  the  exact  values  from  Eq.  (13)  for  the  two  extremes  of  optical  thickness,  0.244  and  12.0, 
in  Fig.  3.  While  the  exponential  approximation  appears  to  give  quite  accurate  results  for 
both  cases,  there  is  no  great  computational  advantage  associated  with  the  approximation. 


3.2  CYLINDRICALLY  SYMMETRIC  CASE 

Analytical  solutions  which  satisfy  the  RTE  in  all  possible  directions  have  also  been  obtained 
for  the  case  of  cylindrically  symmetric  geometry  (Refs.  4 — 5).  These  solutions  represent  the 
radiation  field  inside  an  infinitely  long  circular  cylinder  having  a  black  wall  and  containing 
a  gas  with  temperature  and  absorption  coefficient  profiles  which  depend  on  the  radial 
coordinate  only.  This  configuration  is  a  good  approximation  to  the  AEDC  high-pressure 
arc  heaters.  The  radiant  heat  flux  profile  in  this  configuration  is  given  by  Eq.  (24)  of  Ref, 
4,  which  is  similar  to.  but  more  complex  than  Eq.  (13)  above,  and  contains  exponential  integral 
functions 


which  are  equivalent  to  the  En(Z)  family  (Eq.  1 1)  appearing  in  the  plane-parallel  solution. 

The  arc  heater  radiation  field  was  represented  by  this  model  in  the  Aerotherm  ARCFLOW 
program  (Ref.  5),  except  that  as  in  Eq,  (14),  simple  exponential  functions  were  used  to 
approximate  the  D„(Z)  integral  functions:' 

Di  (Z)  =  —  e-5^4z  and  Dj  (Z)  =  4  (17) 

16  4 

The  McDonnell-Douglas  SWIRLARC  program  (Ref.  6)  contains  the  same  radiation  model. 
Since  it  satisfies  the  RTE  exactly  in  all  directions,  except  for  the  exponential  approximation, 
Eq.  (17),  it  is  considered  to  be  the  most  accurate  procedure  available  for  calculation  of  radiation 
heat  transfer  in  arc  heaters.  It  provides  a  standard  of  comparison  for  approximate  solutions. 
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4.0  APPROXIMATE  SOLUTIONS 

The  complexity  of  the  integro-differential  equation  system,  Eqs.  (8a-b),  has  limited 
analytical  solutions  to  cases  in  which  there  is  a  single  independent  variable  of  position.  Such 
solutions,  however,  cannot  represent  the  radiation  field  associated  with  nonaxisymmetric 
temperature  distribution  such  as  frequently  present  in  a  high-pressure  arc  heater  when  the 
electric  discharge  is  displaced  from  the  constrictor  centerline.  For  this  reason,  approximate 
solutions  which  do  not  satisfy  the  RTE  in  all  directions  of  propagation  have  been  investigated. 

4.1  DIFFUSION  APPROXIMATION 

In  one  approximate  method,  the  behavior  of  the  radiation  Held  is  represented  by  a  system 
of  (only)  partial  differential  equations  of  radiant  heat  flux  and  the  space-integrated  intensity: 

V  .  qR  =  -  K  [lo  -  4<y  T<]  (18) 

V  lo  =  -  3K  qR  (19) 

where 

lo  =  jj”  I  (f.  0  dO  (20) 

This  approximation  of  the  radiation  Held  by  simple  partial  differential  equations  is  referred 
to  as  the  **differential**  or  “diffusion**  approximation  (since  Eq.  (19)  has  the  same  form 
as  the  molecular  diffusion  equation.  Appendix  A.  Sec.  3.0).  These  equations  may  be  developed 
in  several  different  ways,  but  are  shown  in  Appendix  A  to  be  most  appropriately  the  first 
approximation  to  the  radiant  field  expanded  as  an  inHnite  series  of  spherical  harmonic 
functions.  This  series  can  be  made  arbitrarily  close  to  an  exact  solution  for  any  radiation 
field  (Ref.  7). 

Extensive  application  of  this  approximation  in  the  field  of  neutron  transport  theory  has 
demonstrated  that  an  appropriate  (and  approximate)  boundary  condition  is: 

[lo  +  mqR]^„  =  4a  T*^  (21) 

where  the  constant  m  may  have  values  of  y/f  or  2.0  (Appendix  A,  Sec.  2.0). 

4.2  COMPARISON  WITH  EXACT  ANALYTICAL  SOLUTIONS 

The  following  comparisons  were  made  with  the  gas  temperature  specified,  i.e.,  no  coupling 
between  the  gas  energy  equation  and  the  RTE. 


10 


AeOC-TR-91-10 


4.2.1  Plane-Paniild  Walls 

The  radiation  field  between  two  plane-parallel  black  walls  obtained  by  use  of  the  diffusion 
approximation,  Eqs.  (18-19)  (partial  differential  equations  only),  is  compared  with  an  exact 
solution  obtained  by  solution  of  the  integro-differential  system,  Eq.  (IS),  in  Fig.  4.  The  radiant 
heat  flux  distribution  is  given  by  the  same  gas  field  conditions  as  in  Fig.  3,  a  linear  temperature 
distribution  from  25,000  to  5,000°R  and  an  absorption  coefficient  of  24.4  ft~ The  solutions 
are  shown  for  two  extremes  of  optical  thickness,  0.214  and  12.2.  The  diffusion  approximation 
agrees  well  with  the  exact  result,  especially  for  the  optically  thick  case.  Two  values  of  the 
boundary  condition  parameter,  m,  Eq.  (21),  were  used,  the  Mark  value  of  \/3  and  the  Marshak 
value  of  2.0  (Ref.  7).  The  Marshak  value  gives  the  better  agreement. 

4.2.2  Axlsymmetric 

The  radiation  field  for  an  axisymmetric  geometry  approximately  the  scale  of  the  AEDC 
HEAT-Hl  arc  heater  is  shown  in  Fig.  5.  A  linear  gas  temperature  varying  from  25,000'’R 
at  the  centerline  to  5,000°R  at  the  wall  was  used  as  was  a  constant  absorption  coefficient 
of  24.4  ft~L  Exact  solutions  of  Kesten  (Ref.  4),  Chiba  (Ref.  8  ),  and  Nicolet  (Ref.  S)  are 
compared  with  the  present  ‘‘diffusion  approximation”  calculation.  Agreement  was  found 
to  be  good,  especially  at  the  wall. 

5.0  DIFFUSION  APPROXIMATION  IN  TWO  DIMENSIONS 

To  demonstrate  the  utility  of  the  diffusion  approximation  for  multidimensional  fields 
hot  amenable  to  exact  analytical  solutions,  a  two-dimensional  field  in  a  circular  cross  section, 
somewhat  representative  of  the  AEDC  HEAT-Hl  arc  heater  was  solved  using  Eqs.  (18)  and 
(19)  and  the  Marshak  value  of  m  =  2.0  in  the  boundary  condition,  Eq.  (21).  A  gas  temperature 
distribution  representative  of  an  arc  displaced  from  the  centerline  was  assumed:  25,000°R 
peak  temperature  located  at  40  percent  of  the  radius  of  a  2.4-ln.-diam  cross  section  with 
linear  variation  from  this  peak  to  a  uniform  wall  temperature  of  S,000°R.  The  temperature 
distribution  was  thus  decoupled  from  the  radiation  equation.  The  absorption  coefficient  was 
assumed  to  have  a  constant  value  of  24.4  ft~^  independent  of  the  temperature,  and  the  wall 
was  considered  black. 

5.1  SOLUTION  PROCEDURE 

A  Cartesian  coordinate  system  (x,  y)  was  used  to  define  the  circular  cross  section  with 
a  uniformly  spaced  20  x  20  grid.  Variations  in  the  streamwise  direction  (z  coordinate)  were 
neglected.  By  eliminating  the  heat  flux  from  Eqs.  (18-19),  a  single  partial  differential  equation 
was  attained  for  the  space-integrated  intensity  lo¬ 
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The  nuinerical  solution  of  this  equation  for  Ig  was  carried  out  in  an  iterative  manner  as 
follows: 


a.  The  y  derivative  in  Eq.  (22)  is  evaluated  from  the  previous  iteration;  hence,  it 
is  grouped  with  the  source  term  (to  start  the  process  the  initial  iteration  is  taken 
as  Iq  =  constant).  Equation  (22)  then  becomes  only  a  function  of  x  and  is  solved 
with  an  implicit  finite  difference  routine  (Ref.  9)  along  each  grid  line  of  y  = 
constant. 

b.  The  coordinate  system  is  then  rotated  90  deg,  effectively  interchanging  the  roles 
of  X  and  y,  and  the  procedure  outlined  above  in  (a.)  is  repeated  in  an  iterative 
manner. 

c.  The  90-deg  rotations  are  continued  until  the  solution  converges. 


The  boundary  condition,  Eq.  (21),  is  applied  in  a  direction  normal  to  the  boundary,  since 
it  was  developed  for  planar  fields  in  which  the  heat  flux  is  always  normal  to  the  boundaries. 
By  use  of  Eq.  (19),  it  is  recast  in  terms  of  the  integrated  intensity  I®  only: 


lo 


m 

3K 


(23) 


where  n  is  in  the  direction  of  the  inward  normal  to  the  circular  boundary.  In  the  numerical 
solution,  this  boundary  condition  is  applied  at  non-nodal  points  on  the  boundary  by 
interpolation  between  grid  points.  The  final  step  in  the  solution  is  to  obtain  the  radiant  heat 
flux  distribution  from  the  gradient  of  lo,  Eq.  (19). 

5.2  NUMERICAL  RESULTS 

For  an  initial  demonstration  of  the  method  in  two  dimensions,  the  symmetric  temperature 
distribution  on  which  Fig.  5  was  based  was  used  to  calculate  a  radiation  distribution  in  a 
2.4-in.-diam  circle.  The  two-dimensional  calculation  duplicated  the  results  of  the  one¬ 
dimensional  calculation  shown  in  Fig.  S.  Hence,  as  previously  noted,  the  diffusion  results 
agree  well  with  those  of  Kesten  (Ref.  4),  and  others  cited  on  Fig.  5. 

The  same  procedure  was  then  used  to  obtain  the  radiation  field  for  the  nonsymmetric 
temperature  distribution  with  the  peak  temperature  displaced  a  distance  equal  to  40  percent 
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of  the  radius  in  the  positive  y-direction.  In  this  first  approximation  the  gas  temperature  was 
assumed  to  vary  linearly  between  the  maximum  value  and  the  wall  value.  The  results  of  these 
calculations  are  shown  in  Figs.  6a-c.  The  radiant  heat  flux  shown  along  the  Y-axis  in  Fig. 
6a  strongly  reflects  the  effects  of  the  asymmetric  gas  temperature.  Most  notable  is  that  the 
radiant  flux  to  the  wall  nearest  the  maximum  gas  temperature  (i.e.,  arc  centerline)  is  over 
four  times  higher  than  on  the  far  wall.  It  will  be  shown  later  that  the  assumption  of  constant 
absorptivity  greatly  overestimates  the  magnitude  of  this  effect.  The  radiant  flux  along  and 
in  the  direction  of  the  X  axis,  (Fig.  6b),  is  symmetric  about  the  Y-axis  (because  the  gas 
temperature  is  symmetric);  however,  the  maximum  flux  occurs  at  a  larger  radius  than  it  does 
in  the  axisymmetric  case.  The  radiant  flux  along  a  4S-deg  line  (X  =  Y)  (Fig.  6c)  simply  shows 
a  blending  of  the  two  extremes  (Figs.  6a-b).  The  radiant  heat  flux  at  the  wall  is  depicted  in 
Fig.  6d  at  various  positions  around  the  circular  cross  section.  The  relationship  between  the 
location  of  the  temperature  peak  and  the  maximum  radiant  wall  heat  flux  is  readQy  apparent. 

The  solution  for  the  two-dimensional  demonstration  problem  provides  a  high  degree  of 
confidence  that  the  “diffusion  approximation”  can  be  used  in  a  three-dimensional  arc  heater 
code  currently  being  developed  since  the  primary  gradients  will  be  in  the  cross-sectional  plane. 

6.0  DIFFUSION  APPROXIMATION  APPLIED  TO  SWIRLARC  CODE 

To  gain  some  experience  in  solving  the  gas  energy  equation  and  the  “diffusion 
approximation”  RTE  in  their  respective  coupled  forms,  it  was  decided  to  incorporate  the 
“diffusion  approximation”  into  the  axisymmetric  SWIRLARC  code  (Ref.  6). 

6.1  SOLUTION  PROCEDURE 

The  approximate  radiation  model  represented  by  Eqs.  (18),  (19),  and  (21)  was  added  to 
the  system  of  equations  in  the  SWIRLARC  code,  replacing  the  relations  representing  the 
exact  analytical  solution  of  the  radiative  transfer  equation.  The  numerical  solution  procedure 
used  in  SWIRLARC  is  an  explicit  marching  procedure.  To  ensure  that  the  solution  procedure 
was  adequate  for  the  RTE,  an  implicit  solution  procedure  with  one-half  the  radial  mesh  size 
was  used  for  it.  The  respective  coupling  terms  between  the  gas  energy  equation  and  the  RTE 
were  included.  The  thermal  property  curve  fits  which  are  built  into  SWIRLARC  (including 
those  for  absorption  coefficient)  were  used  for  all  calculations. 

6.2  NUMERICAL  RESULTS 

A  pair  of  solutions  was  first  obtained  using  the  standard  SWIRLARC  code  with  both 
13  and  25  radial  mesh  points  to  determine  the  influence  of  mesh  size  on  details  of  the  solution, 
notwithstanding  the  fact  that  previous  calculations  had  shown  that  overall  results  are  not 
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sensitive  to  the  radial  mesh  size.  These  two  calculations  were  made  for  the  HEAT-Hl  geometry 
at  conditions  of  P*  =  22  atm,  =  3,570  Btu/lbm,  I  =  1,200  amp,  and  Z  =  1.2  m  in 
a  S.08-cm  (2-in.)-diam  bore.  The  calculated  radiant  heat  flux  and  temperature  distributions 
are  shown  on  Fig.  7a.  Although  the  maximum  radiant  heat  flux  inside  the  radiant  field  is 
considerably  higher  for  the  13-point  mesh,  the  values  at  the  wall  agree  within  7  percent.  This 
does  not  seem  to  be  a  significant  difference  in  light  of  the  uncertainties  in  the  absorption 
coefficient  for  high-temperature  air  to  be  discussed  in  Sec.  8.0.  The  effect  of  mesh  size  on 
the  radial  temperature  distribution  is  almost  negligible,  except  near  the  centerline  of  the  flow 
field. 

Solutions  with  the  diffusion  approximation  of  radiation  in  the  SWIRL  ARC  code  were 
obtained  with  a  13-point  radial  mesh  for  solving  the  gas  dynamic  equations  and  a  25-point 
mesh  for  the  radiation  equations.  The  resulting  calculated  radiant  heat  flux  distribution  for 
the  same  Po  =  22  atm  case  (Fig.  7a)  was  found  to  fall  between  the  13-point  and  2S-point 
calculation  of  the  standard  SWIRLARC  code,  and  the  radial  temperature  distributions  were 
essentially  unchanged. 

Similar  calculations  (Fig.  7b)  were  also  made  for  the  case  of  ?„  =  115  atm.  Ho  =  3,200 
Btu/lbm,  I  =  1,600  amp,  and  Z  =  1.6  m.  For  this  high-pressure  condition,  differences  in 
calculated  radiant  heat  flux  and  temperature  distributions  introduced  by  the  diffusion 
approximation  were  even  less  than  for  the  lower  pressure  case. 

It  appears  that  at  least  for  axisymmetric  flow  and  radiation  fields,  the  diffusion 
approximation  of  radiative  energy  transport  can  be  readily  coupled  to  a  gas  dynamic  calculation 
and  will  produce  temperature  and  radiant  energy  transfer  distributions  that  are  reasonably 
close  to  those  obtained  with  an  exact  model  of  the  geometries  of  the  radiation  field. 

7.0  INTERPRETATION  OF  EXPERIMENTAL  DATA 

During  a  test  in  the  AEDC  HEAT-Hl  arc  heater  during  FY90,  high-speed  motion  pictures 
were  made  through  a  quartz  side  window  in  the  constrictor  wall.  Movies  were  obtained  of 
the  image  of  a  short  segment  of  the  arc  (defined  by  a  narrow  slit)  as  it  moved  relative  to 
the  centerline  reference.  A  schematic  of  the  HEAT-Hl  arc  heater  showing  the  general  heater 
configuration  and  the  locations  where  the  window  segment  was  installed  is  presented  in  Fig. 
8.  Data  were  obtained  at  50  and  80  atm  with  the  window  located  in  Modules  2,  5,  and  8. 
Short  sections  of  movie  film  (16  msec)  taken  near  the  downstream  end  of  the  heater  (i.e.. 
Module  2,  where  heat-transfer  rates  are  the  highest)  are  shown  in  Fig.  9.  Excursions  of  the 
arc  image  from  the  centerline  reference  are  presumably  the  result  of  (nearly)  circular  motion 
of  the  arc  column  inside  the  heater,  and  the  width  of  the  trace  is  an  approximate  measure 
of  the  diameter  of  the  arc  column  cross  section. 
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7.1  MEASURED  VERSUS  .CALCULATED  HEAT  FLUX 

Heat  flux  to  the  arc  heater  constrictor  wati  is  plotted  versus  axial  position  in  Fig.  10a 
for  Po  =  52,4  atm,  Ho  =  2,700  Btu/lb,  and  in  Fig.  10b  for  Po  =  82.5  atm,  Ho  =  3,120 
Btu/lb.  The  experimental  data  are  shown  as  bars  which  represent  the  average  heat  flux  for 
each  of  the  eight  modules,  based  on  measured  temperature  rise  in  cooling  water.  The  predicted 
curves  shown  were  calculated  with  the  SWIRLARC  code.  Note  that  the  predictions  for 
Modules  7  and  8  are  suspect  because  the  arbitrary  initial  profiles  used  to  start  the  SWIRLARC 
code  have  not  yet  “washed  out,“  and  the  drop  in  heat  flux  between  1 .5  and  1 .6  m  is  caused 
by  the  addition  of  the  final  increment  of  air  which  in  reality  is  added  in  the  electrode  region 
downstream  of  the  constrictor  (not  modeled  in  SWIRLARC). 

The  heat  load  in  the  upstream  end  of  the  heater  is  primarily  radiative.  The  radiative  load 
is  nearly  constant  while  the  convective  load  increases  in  the  downstream  direction  as  the  gas 
velocity  increases.  As  a  check,  the  radiative  component  was  calculated  each  of  the  following 
three  ways  with  the  same  result; 

1 .  SWIRLARC  code  with  the  integrated  radiation  model; 

2.  SWIRLARC  code  with  the  differential  radiation  model,  axially  symmetric;  and 

3.  Differential,  2-D,  uncoupled  radiation  code  using  a  symmetric  gas  temperature 
profile  from  SWIRLARC. 

Overall,  there  is  reasonably  good  agreement  between  the  measured  and  predicted,  wall 
heat  flux. 

7.2  EFFECT  OF  AN  OFF-CENTER  ARC  ON  THE  WALL  HEAT  TRANSFER 

In  the  previous  section  the  two-dimensional  diffusion  radiation  code  was  shown  to  give 
results  consistent  with  other  predictive  methods  for  symmetrical  fields  and  experimental  data 
averaged  over  each  module.  The  radiation  code  was  next  used  to  investigate  the  effect  of 
a  fluctuating  arc  location  (shown  in  Fig.  9)  on  the  local  wall  heat  transfer.  The  following 
modifications  were  made  to  the  2-D  radiation  code. 

a.  The  grid  was  increased  to  (80  X  80). 

b.  The  pressure  and  temperature  dependence  of  the  absorption  coefflcient,  K,  was 
included. 

c.  The  basic,  on-centerline  arc  temperature  profile  was  taken  from  SWIRLARC  and 
approximated  by  four  piecewise  linear  segments. 
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d.  The  off-centerline  arc  temperature  profile  was  obtained  by  shifting  the  break  points 
between  the  linear  segments  by  a  proportional  amount.  This  process  is  illustrated 
in  Fig.  11  for  two  linear  segments. 

Radiant  heat  flux  calculations  were  carried  out  for  arc  displacements  up  to  90  percent 
of  the  constrictor  radius  with  reservoir  conditions  representative  of  HI  operation  at  SO  and 
80  atm  and  also  for  a  projected  operation  at  190  atm.  The  radiant  flux  profiles  were 
qualitatively  similar  to  those  shown  in  Fig.  6.  The  ratios  of  maximum  and  minimum  wall 
heat  flux  (near  and  far  wall)  to  the  symmetric  arc  value  are  shown  in  Fig.  12  as  a  function 
of  the  arc  displacement  for  the  reservoir  conditions  noted.  The  assumption  of  constant 
absorptivity  used  in  the  textbook  problem  of  Fig.  6  (K  =  24.4  ffi,  for  100  atm,  25,000®R) 
greatly  exaggerates  the  effect  of  the  arc  displacement.  However,  linear  fits  to  the  calculated 
data  indicate  a  50-percent  increase  in  the  near-wall  radiant  heat  flux  as  the  arc  approaches 
the  wall  at  195  atm.  The  decrease  in  the  far-wall  radiant  heat  flux  is  of  the  same  order  as 
the  increase  at  the  near  wall. 

7.3  WALL  HEAT-TRANSFER  ANALYSIS  WITH  UNSTEADY  BOUNDARY 

CONDITIONS 

A  one-dimensional,  unsteady  numerical  model  for  conduction  heat  transfer  through  the 
0.080-in.  arc  heater  constrictor  wall  was  developed.  The  water-side  wall  temperature  was 
taken  as  a  constant  650“F.  This  is  approximately  100®F  above  the  saturation  temperature 
at  1 ,000  psia  (cooling  system  pressure)  and  represents  a  practical  limit  on  the  water-side  wall 
temperature  to  avoid  burnout.  This  assumption  will  eventually  be  improved  with  the  results 
of  a  boiling  heat-transfer  experiment  currently  underway  at  AEDC.  The  unsteady  air-side 
boundary  condition  is  formulated  by  first  assuming  that  the  arc  location  (Fig.  9)  can  be 
represented  by  a  linear  oscillation,  a  fairly  poor  assumption  since  the  path  is  probably  three- 
dimensional.  For  simplicity,  the  arc  position  is  approximated  by  a  sine  wave  with  period 
and  amplitude  matched  to  the  film  record  (Fig.  9).  The  instantaneous  location  of  the  arc 
is  then  used  with  the  results  shown  in  Fig.  12  to  determine  the  radiant  wall  heat  transfer. 
Finally,  the  convective  component  is  added  to  the  radiative  component  to  obtain  the  total 
heat  transfer  which  serves  as  the  air-side  boundary  condition  for  the  conduction  calculation. 

Results  of  a  calculation  for  a  case  of  SO  atm  and  3,120  Btu/lbm  are  shown  in  Fig.  13a 
for  a  period  of  0.001  sec  and  an  amplitude  of  90  percent  of  the  constrictor  radius.  This  is 
a  typical  value  of  the  period  (Fig.  9)  and  a  maximum  value  of  the  amplitude.  The  frequency 
is  sufficiently  high  so  that  the  water-side  heat  flux  is  not  affected  and  the  air-side  wall 
temperature  rises  only  a  few  degrees.  The  possibility  of  the  arc  remaining  near  the  wall  for 
an  extended  period  of  time  was  investigated  by  using  a  period  of  0.1  sec.  The  results  shown 
in  Fig.  13b  show  significant  rises  in  both  the  air-side  wail  temperature  and  the  water-side 
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heat-transfer  rate.  However,  the  maximum  heat  flux  of  4,000  Btu/ft^  sec  at  80  atm  is 
sufficiently  low  that  wall  failure  is  unlikely,  even  if  the  arc  remains  close  to  the  wall  for  an 
indefinite  period. 

A  case  at  195  atm,  3,520  Btu/lbm  was  also  investigated.  A  period  of  0.001  sec  (taken 
from  the  50-  and  80-atm  case)  and  an  amplitude  of  90  percent  of  the  constrictor  radius  was 
used  to  obtain  the  results  shown  in  Fig.  14a.  Again,  the  frequency  of  the  arc  fluctuation 
is  sufficiently  high  so  that  the  water-side  of  the  constrictor  wall  does  not  respond.  The  air- 
side  wall  temperature  responds  only  slightly,  rising  about  30°F.  The  case  of  the  arc  remaining 
near  the  wall  for  an  extended  period  was  also  investigated  using  a  period  of  0.1  sec.  The 
results  shown  in  Fig.  14b  show  significant  rises  in  both  the  water-side  heat  transfer  and  the 
air-side  wall  temperature.  The  wall  temperature  rises  to  1,830®F,  which  is  still  below  the  melting 
point.  The  water-side  heat-transfer  rate  peaks  at  10,600  Btu/ft^  sec,  which  is  near  the  value 
at  which  burnout  is  thought  to  occur.  A  response  time  of  sorts  may  be  obtained  by  looking 
at  the  phase  shift  in  Fig.  14a.  For  the  air-side  wall  temperature  this  is  0.0088  sec,  and  0.015 
sec  for  the  water-side  heat  transfer. 

Finally,  the  temperature  profiles  through  the  copper  wall  at  the  maximum  wall  temperature 
are  shown  in  Fig.  15.  If  the  arc  remains  on  the  geometric  centerline,  the  profile  is  linear. 
For  the  low-frequency  fluctuation  (X  =  0.1  sec)  the  deviation  from  linearity  is  imperceptible, 
while  for  the  high-frequency  fluctuation  (X  =  0.001  sec)  the  nonlinearity  begins  to  show  up. 
The  tensile  strength  of  copper  for  the  X  =  0.1  sec  case  is  also  plotted  in  Fig.  15.  Although 
melting  has  not  occurred,  the  strength  of  the  wall  on  the  air  side  is  quite  low. 

8.0  RADIATION  PROPERTIES  OF  HIGH-TEMPERATURE  AIR 
8.1  TWO-BAND  MODEL 

An  attempt  was  made  to  resolve  the  differences  in  presentation  of  radiative  properties 
of  high-temperature  air  in  a  number  of  references.  In  Ref.  5,  air  was  characterized  as  having 
two  gray  emission  bands,  a  low-energy  band  with  0  - 10.5  ev  photon  energy  and  a  high-energy 
band  with  10.5  -  100  ev  photon  energy. 

Between  2,000  and  30,000*K  the  low-band  absorption  coefficients  are  Kl  =  10“ ^  to  4 
X  10“’  cm~’  and  the  high-band  absorption  coefficients  are  Kh  =  1.0  to  1  x  10^  cm-’. 
Although  the  absorption  coefficients  are  much  larger  for  the  high-energy  baiid,  the  fraction 
of  total  blackbody  intensity  is  also  very  small  for  this  band,  <  1  percent  for  temperatures 
up  to  11,000°K.  It  is  therefore  possible,  under  certain  conditions,  to  obtain  a  good 
approximation  to  total  radiant  intensity  using  only  the  low-band  absorption  coefficient,  K^, 
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lit®  -  Bt  (1  -  ^  (1  -  (24) 


The  error  made  in  use  of  this  relation  may  be  assessed  by  evaluating  the  ratio: 
Itotal  _  II  +  Ih  _  Bl(1  -  Bh  (l  - 

IiND  IiND  Bt  (]  — 


=  Wl  +  Wh 


(I  -  e-’^HX) 
(1  -  e-'^L’^) 


1  +  Wh 


(I  -  e-^HX) 

(1  - 


(25) 


where  Bl,  Bh  are  blackbody  intensities  for  the  two  bands  and  Wl,  Wh  are  blackbody  band- 
weight  factors. 


A  plot  of  this  ratio  is  given  in  Fig.  16  for  pressures  of  1,  10,  60,  and  100  atm,  at  1-cm 
depth  of  gas  and  for  O.S  and  2  cm  at  10  atm.  The  emission  from  air  is  clearly  dominated 
by  the  low-energy  band  for  sufficiently  high  pressures  and  temperatures  less  than  15,000  K. 
The  low-band  absorption  coefficient  may  be  used  with  good  accuracy  at  T  <  12,000  K,  P 
a  60  atm,  and  X  ^  1  cm,  all  of  which  conditions  prevail  in  the  high-pressure  AEDC  arc 
heaters  under  normal  operating  conditions.  For  pressures  of  the  order  of  10  atm  and  below, 
however,  the  errors  can  become  quite  large,  &  100  percent. 

8.2  ABSORPTION  COEFFICIENT  DATA 


Radiation  properties  of  high-temperature  air  are  available  in  a  large  number  of  references. 
A  comparison  is  complicated  because  of  the  diverse  forms  in  which  the  data  are  presented 
and  also  because  the  thermodynamic  states  are  generally  not  matched.  Integrated  hemispheric 
emission  coefficients  are  given  in  Refs.  10  -  12,  spectral  and  integrated  steradiance  in  Refs. 
13  -  14,  and  integrated  Planck  absorption  coefficients  in  Refs.  5  and  11.  Most  of  these 
references  rely  on  data  from  precedent  original  references.  In  some  cases,  the  emission  data 
are  broken  down  by  individual  radiating  species  in  the  air.  For  comparison  of  these  data, 
all  were  converted  to  total  (integrated)  absorption  coefficient,  interpolated  in  some  cases  to 
match  specific  thermodynamic  states.  The  absorption  coefficient  at  I -atm  density  and 
temperature  from  2,000  to  13,000  K  is  given  in  Fig.  17a,  from  data  obtained  in  the  cited 
references.  On  Fig.  17b,  the  absorption  coefficient  is  given  as  a  function  of  density  at  T  = 
6,000  K.  The  low-band  coefficients  of  Nicolet  {Ref.  5),  which  are  also  used  in  the  SWIRLARC 
program,  (Ref.  6),  agree  with  data  from  other  sources  within  a  factor  of  2.0  to  3.0. 
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9.0  CONCLUSIONS 

Computational  exercises  have  been  used  to  determine  whether  or  not  a  simplified  radiation 
model  (diffusion  approximation)  can  be  used  for  the  calculation  of  high-temperature  gas 
flows  in  which  the  radiation  fields  are  not  sufficiently  symmetric  to  be  amenable  to  ^act 
analytical  solution  of  the  radiative  transfer  equation.  Comparisons  of  radiation  fields  calculated 
exactly  and  with  various  approximations  were  made  for  simple  configurations  for  which 
analytical  solutions  are  available.  Also,  simple  axially  symmetric  flows  of  hot  gases  were 
successfully  calculated  using  approximate  radiation  field  equations  coupled  to  the  gas  dynamic 
conservation  equations. 

1.  Based  on  these  results,  there  appears  to  be  a  reasonable  basis  for  applying  the 
diffusion  approximation  of  radiation  in  flow  calculations  in  which  the 
temperature  and  radiation  fields  are  nonsymmetrical,  such  as  in  an  electric  arc 
heater  with  the  discharge  displaced  from  the  centerline. 

2.  Total  wail  heat-transfer  rates  measured  in  the  constrictor  of  the  AEDC-H 1  arc 
heater  at  50  -  80  atm  pressure  and  2,700  -  3,100  Btu/lbm  enthalpy  agreed  quite 
well  with  predictions  made  by  the  SWIRLARC  computer  program.  The  radiant 
flux  dominates  at  the  upstream  end  of  the  constrictor. 

3.  An  unsteady  heat-transfer  analysis  of  the  arc  heater  constrictor  walls  revealed 
that  the  observed  fluctuations  in  arc  position  were  of  sufficiently  high  frequency 
that  the  walls  will  not  respond  to  the  fluctuating  heat  input. 
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Observation  Direction  of 


a.  Radiant  Intensity  in  an  emitting  and  absorbing  gas 


Optical  Thickness,  ii  (x)  a 

o 

b.  One-dimensional  case  (radiant  field  between  two  parallel,  infinite  black  walls) 
Figure  1.  Geometry  considerations  for  radiant  energy  transfer  in  an  emitting  and 
absorbing  gas. 
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Current 

Kesten  (Ref.  4)  - 
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Thickness 
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A  E  0.60  in. 

1.22 

A  E  1.20  in. 

2.44 

A  s  3.00  in. 

6.10 

A  s  6.00  in. 

12.20 

Exact  Solution  of  Radiative  Transfer  Equation,  Linear  T(X),  25,000“  -  5,000“R, 
Black  Walls,  Absorption  Coefficient  =  K  s  24.1  ft-i . 


0  0.20  0.40  0.60  0.80  1.00 

Figure  2.  Radiation  heat  flux  between  plane  walls  —  optical  thickness  effects. 
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Figure  3.  Rsdiation  heat  flux  between  plane  walls  —  exponential  approximation 
compared  with  exact  solution. 
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Figure  4.  Radiation  heat  flux  between  plane  walls  —  diffusion  approximation  compared 
with  exact  solution,  linear  T(X),  25,000  —  5,000°R,  black  walls,  absorption  coef¬ 
ficient  =  24.4  ft”*. 
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Figure  5.  Radiation  flux  in  long  circular  cylinder  —  comparison  of  diffusion  approximation 
with  exact  solutions. 
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r/R 

c.  qjfy  along  X  =  Y  diagonal 

Figure  6.  Radiation  heat  flux  distributions  In  long  circular  cylinder  with  asymmetric 
temperature  distribution,  temperature  peak  at  r/R  »  0.40,  R  =  1.2  in. 
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Note:  All  Dimensions  in  Inches 


Is! 


16.200  Btu/ft2sec 


d.  Distribution  of  mdiant  beat  flux  at  with  wall  asymmetric  temperature  distribution: 
ZStOOO^R  at  r/R  «  0.4,  T^aii  =  5,000'’R,  linear,  absorption  coefficient  =  24.4  ft"  ^ 

Figure  6.  Concluded. 
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a.  Po  =  22  atm.  Ho  =  3,570  Btu/lbm,  I  =  1,200  amps,  Z  «  1.2  m,  d  =  5.08  cm 
Figure  7.  Temperature  and  radiant  heat  flux  distribution  in  SWIRLARC  solutions. 
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n,  I  =  1,600  amps,  Z  =  1.6  m,  d  =  5.08  cm 
7.  Concluded. 


oi-ie-ui-oaav 


( 


I 


Anode  “ 
Section 


nodules  8-1- 


Cathode 

Section 


Figure  8.  Window  segment  locations  —  HEAT-Hl  arc  heater. 
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DConstrittor 

=  2.0  in. 


Centerline  Reference  Mark  (Typ.) 


Chamber  Pressure:  52  atm 

Timing  Marks  Enthalpy:  2,730  Btu/lb 

(One  every  1  msec) 


Note:  The  camera  exposure  setting  was  changed  for  the 

82-atm  run  to  prevent  overexposure,  resulting  in 
a  difference  in  arc  intensities. 


Dconstrictor 

=  2.0  in. 


Chamber  Pressure:  82  atm 
Enthalpy:  3,120  Btu/lb 


Figure  9.  Arc  visualization  at  two  different  chamber  pressures. 
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qwall,MW/m2 


Distance  from  the  Upstream  End  of  the  Arc  Heater,  Z,  m 
a.  52.5-atiii  runs 


Figure  10.  Heat  flux  measurements  in  the  HEAT-Hl  arc  heater  and  comparison  with 
SWIRLARC  predictions. 
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Po  s  82.5  atm.  Hq  »  3.120  Btu/lbm 
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Distance  From  the  Upstream  End  of  the  Arc  Heater,  Z.  m 
b.  82.5-atm  ran 
Figure  10.  Concluded. 
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Figure  11.  Postulated  gas  temperature  distribution  resulting  from  an  off-centerline  arc  location. 
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Temperature  Profiles  Based  on  SWIRL  ARC  Solutions  for: 

•  Ho  =  3,S20  Btu/lb,  pQ  e  195  atm.  Z  =  1.4  m 
A  Ho  -  3,120  Btu/lb.  Po  =  82.5  atm,  (HI >078-01 0),Z  =  1.4m 
O  Ho  "  2,730 Btu/lb, Po  »  52.4 atm, (HI -078-005), Z  =  1.4m 
□  Ho  c  2,740  Btu/lb,  Po  »  52.4  atm,  (HI -078-009),  Z  -  0.25  m 


/ 

/ 


Arc  Displacement  from  Centerline,  percent  of  Radius 
Figure  12.  Radiant  wall  heat  flux  modification  due  to  off-centerline  arc  location. 
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Air  Side  Water  Side 


Figure  15.  Predicted  temperature  and  strength  profiles  in  constrictor  wall  at  maximum 
wall  temperature. 
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I|ND  =  Bt(1  ■  e-KLX) ;  Bj  *  aT^/n ,  Ku  a  Abs.  Coefficient.  10.5  ev  Band 


Figure  16.  Ratio  of  total  intensity  to  indicated  intensity  versus  temperature  (air). 
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B.  Zel  Dovich/Keck/Meyerott,  Ref.  14 

C.  Martin/Breene.Nardone,  Ref  13 

D.  Renner/Gilmore.  Ref.  1 1 


2.000  4.000  6,000  iB,000  10,000  12,000  14,000 

Temperature,  K 


a.  q/ Po  —  1.0  (EphoiM  ”  0  “  10.5  ev) 

Figure  17.  Absorption  coefficient  of  air  versus  temperature  and  density. 
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b.  T  =  6,000  K,  Ephotoo  =  0  -  10.5  ev 
Figure  17.  Concluded. 
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APPENDIX  A 

SOLUTION  OF  RADIANT  ENERGY  HELD  BY  SPHERICAL  HARMONIC  FUNCTIONS 

1.0  INFINITE  SERIES  EXPANSION 


As  described  in  Ref.  7,  any  arbitrary  radiation  field  may  be  exactly  specified  by  an  infinite 
series  of  spherical  harmonic  functions.  In  three  dimensions,  the  radiative  transfer  equation 
(Eq.  3)  is: 


I  .  VI 


(A-1) 


where  F  is  a  unit  vector  in  the  direction  of  propagation  and  1  is  the  frequency-integrated 
radiant  intensity,  Eq.  (2).  It  is  assumed  that  a  general  solution  to  Eq.  (A-1)  can  be  represented 
by  an  infinite  series: 

00  +f 

•  ^  ,  ?o  n,  Y^IF)  (A-2) 


where  Y^{1)  are  the  spherical  harmonic  functions  of  direction  and  A)°(r)  are  spatially 
dependent  coefficients  to  be  determined  by  boundary  conditions  on  the  radiation  field. 


By  a  standard  procedure  outlined  in  Ref.  7,  the  series  Eq.  (A-2)  is  substituted  in  the  RTE, 
Eq.  (A-1),  followed  by  multiplication  by  the  complex  conjugate  of  y”  integration  over  47r 
solid  angle,  and  imposition  of  orthogonality  and  recurrence  relations.  A  set  of  partial 
differential  equations  is  obtained  having  the  form: 


6 

I 


a  =  1 


fa  (/.m)  ga 
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II1±  \ 
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3*1.3 


3A”±  I  ^ 
3X2 
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-K 


v 


(A-3) 


There  is  one  such  equation  for  each  pair  of  the  I  and  m  indices.  The  fi's  are  Kronecker 
deltas  and  the  Xj  are  Cartesian  coordinates  related  to  r,  $,  spherical  coordinates  by  the 
direction  cosines  of  the  unit  propagation  vector,  t : 


fi  =  X|/r  =  sin3  sin)^;  (3  =  *2/r  =  cos3;  ^3  =  xs/r  =  sin3  cos0  (A-4) 


In  the  limit,  as  I  this  infinite  set  of  simultaneous  partial  differential  equations  is 
completely  equivalent  to  the  integro-differential  equations  of  the  exact  solutions.  An 
approximation  to  this  system  is  obtained  by  terminating  the  series  representation  for  the 
intensity  function  I  after  the  i  =  term.  There  will  always  be  (N  +  1)^  equations  of  the 
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intensity  function  I  after  the  I  =  Nth  term.  There  will  always  be  (N  +  1)^  equations  of  the 

type  of  Eq.  (A-3)  involving  the  same  number  of  A^coefficients.  As  N  increases,  the  system 
of  equations  becomes  a  closer  and  closer  approximation  to  the  infinite  set.  It  is  known  from 
neutron  transport  theory  that  the  first  approximation  (N  =  1)  is  sufficiently  accurate  for 
many  radiation  flelds;  that  is: 

I  (f,  D  -  aJ  yJ  +  A,-'  V-‘  +  aX  +  a;  Y,‘  +  ...  (A-5) 

The  four  equations  of  type  Eq.  (A-3)  in  this  approximation  could  be  solved  directly  for 
the  four  A,  coefHcients,  thus  deHning  the  radiation  field  intensity  for  given  boundary 
conditions.  Instead,  however,  advantage  is  taken  of  the  fact  that  the  space-integrated  intensity 


lo  =  jj'  I  Cr.  0  dQ  (A-6) 

and  the  radiant  flux  vector  q^  are  also  functions  of  the  four  A“coefficients,  so  that  in 
combination  with  the  four  partial  differential  equations  of  type  Eq.  (A-3),  the  following  two 
relations  are  developed: 


3qi  ^  5q2  ^  ^q3 
dxi  dx2  dx3 


-K[U  -  4^1^] 


aio 

dX| 


+ 


d\2 


'2  + 


aio  . 


-3K  [q|ii  +  q2  i2  +  Q3  b] 


(A-7) 

(A-8) 


In  vector  notation,  finally. 


V  .  qR  =  -K  [lo  -  4ffT<]  (A-9) 

VIo=-3KqR  (A-10) 

Equations  (A-9)  and  (A-10)  are  identical  to  the  relations  which  arc  obtained  by  moment 
methods  and  also  referred  to  as  the  diffusion  or  differential  approximations  (Sec.  4.1). 

If  the  appropriate  spherical  harmonic  functions  (Ref.  7)  are  substituted  in  Eq.  (A-S),  the 
radiation  intensity  in  this  approximation  can  be  expressed  as 

I  (r.F)  =  (lo  (f)  +  3  qR  .  r)  (A-11) 


for  general  radiation  fields,  and 

I  (f,  O  =  (lo  (f)  +  3  q  costf) 


(A-12) 
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for  simple  planar  radiation  fields  in  which  B  is  the  angle  between  direction  of  propagation 
and  the  heat  flux  vector  q®.  The  intensity  is  thus  the  spatial  average  value  Io/4t  plus  a 
directionally  dependent  contribution  related  to  the  heat  flux. 

2.0  BOUNDARY  CONDITIONS 

The  boundary  condition  associated  with  the  exact  radiative  transfer  Eq.  (A-1)  is  simply 
that  the  outward  radiant  intensity  at  a  black  wall  is  equal  to  the  blackbody  intensity  of  the 
wall  temperature,  uniform  in  all  directions: 

I*  (;,  r)  = 

If  the  RTE  is  properly  integrated  with  respect  to  distance  from  this  wall  value  over  all 
possible  directions,  with  a  known  gas  temperature  distribution,  the  exact  solution  for  I  (f, 
T)  i&  uniquely  determined.  If  the  spherical  harmonic  expansion  of  the  radiation  field  is  used, 
on  the  other  hand,  this  boundary  condition  is  replaced  by  an  infinite  set  of  boundary  conditions 
obtained  by  substitution  of  the  series  expansion  for  intensity  Eq.  (A-2)  in  (A-13).  In  the  limit 
of  f  —  00,  this  set  of  boundary  conditions  combined  with  Eq.  (A-3)  is  entirely  equivalent 
to  the  combination  of  Eqs.  (A-1)  and  (A-13). 

In  any  finite-order  approximation  of  spherical  harmonics,  however,  the  exact  boundary 
condition  cannot  be  met  because  a  truncated  series  for  I,  such  as  Eq.  (A-3),  cannot  represent 
a  directionally  independent  intensity  at  a  wall,  except  for  the  case  of  zero  radiant  heat  flux. 
An  approximation  to  the  boundary  condition,  Eq.  (A-13),  must  be  adopted  and,  moreover, 
many  different  approximations  are  possible.  The  most  frequently  used  boundary  condition 
approximations  arc  those  developed  in  neutron  transport  theory  for  problems  with  planar 
geometry  and  dependence  upon  a  single  space  coordinate  (Ref.  7): 

1 .  Mark’s  boundary  condition  simply  requires  that  the  radiant  intensity,  I  =  l/4ir 
X  (lo  +  3  q  cos  (?),  at  the  wall  be  equal  to  the  blackbody  intensity  only  in  a  single 
discrete  direction  given  by  0-,  =  cos“‘  (I  v^)=  54.7“. 


2.  Marshak’s  boundary  condition  requires  an  intensity  distribution  at  the  wall  which 
produces  an  integrated  radiant  heat  flux  from  the  wall  equal  to  that  of  a  blackbody. 

These  two  approximate  boundary  conditions  then  relate  integrated  intensity  Iq,  heat  flux 
q,  and  wall  temperature  by 


[lo  +  raq],,.i,  =  4ffT^ 


(A-14) 


with  m  =  V3  for  the  Mark  condition  and  m  =  2.0  for  the  Marshak  condition. 
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Reference  7  states  that  extension  of  this  relation  to  more  general  radiation  fields  is  neither 
straightforward  nor  free  of  ambiguities.  On  the  other  hand,  the  three-dimensional  field 
associated  with  a  variable-temperature  plane  wall  has  a  first-order  Marshak  boundary  condition 
which  is  cited  as: 


lo  (Xi,  0,  X3)  +  2q2  (X|,  0,  X3)  =  4(rT^  (Xi,  0,  Xj)  (A-15) 

only  marginally  different  from  Eq.  (A-14).  Based  on  this  similarity,  it  is  assumed  that  boundary 
conditions  for  other  nonplanar  fields  such  as  those  occurring  in  arc  heaters  may  be  obtained 
from  Eq.  (A-i4)  without  serious  error. 

3.0  THE  DIFFUSION  APPROXIMATION 

The  approximate  solution  to  the  problem  of  radiant  transport  of  energy,  Eqs.  (A-9)  and 
(A-10)  with  appropriate  boundary  conditions,  is  often  referred  to  as  the  “diffusion” 
approximation.  This  is  based  on  the  similarity  of  Eq.  (A-10)  to  the  relation  for  molecular 
diffusion  of  neutral  particles: 

„  dn  Lc  dn  ,,, 

Qdiff  =  —  = - 1 - ^ —  (A-16) 

dx  3  dx 

where  qdirr  is  the  molecular  number  flux,  D  is  the  binary  diffusion  coefficient,  L  is  the  mean 
free  path,  c  is  mean  molecular  speed,  and  dn/dx  is  the  concentration  gradient.  The  analogy 
is  made  apparent  by  converting  Eq.  (A-10)  to  a  photon  flux: 

photon  flux  =  *= - ^ —  V  fnphi'c)  =  V  (no)  {A-17) 

hi*  3Khi.  y  ^  f  3K 

where  c  is  the  speed  of  light.  The  absorption  coefficient,  K,  is  thus  identified  as  the  reciprocal 
of  a  mean  free  path  in  absorption. 
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NOMENCLATURE 


a,b 

A,Y 

B 

c 

D 


Di,D2,D3 


Ei,E2iE3 


h 

H 

i 

[ 

I, 

lo 

K 

f 

r 

L 

m 


Constants,  Eq.  (14) 

Coefficients,  Eq.  (A-2) 

Blackbody  intensity 

Mean  molecular  speed  or  speed  of  light 
Binary  diffusion  coefndent 
Exponential  integral  functions,  Eqs.  (1 1)  and  (16) 
Exponential  functions,  Eq.  (17) 

Exponential  functions,  Eq.  (14) 

Planck’s  constant 
Enthalpy 
Unit  vector 

Frequency-integrated  radiant  intensity  or  current 
Spectral  radiant  intensity  function 
Space  integrated  radiant  intensity 
Absorption  coefficient 
Direction  cosine 

Unit  vector  in  direction  of  propagation 
Mean  free  path 
Constant,  Eq.  (A- 14) 
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n  Direction  normal  to  surface  or  particle  density 

P  Pressure 

q  Heat  flux 

r  Position  vector 

r,0,^  Spherical  coordinates 

S  Distance  in  direction  of  propagation 

T  Temperature 

W  Band  weight  factor 

X,Y,Z  Cartesian  coordinates 

5  Kronecker  delta 

7}  Optical  thickness,  Eq.  (9) 

K  Frequency 

a  Stefan-Boltzman  constant  - 

r  Optical  depth,  Eq.  (S) 

n  Solid  angle 

A  Thickness  (distance  between  plane  walls) 

X  Period,  sec 

Subscripts 

B  Boundary 

H  High-energy  band 
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i  Coordinate  indices  1,  2,  3 

L  Left  or  low-energy  band 

o  Standard  temperature/pressure 

p  Photon 

R  Right 

S  Position  along  direction  of  propagation 

W  Wall 

V  Photon  frequency 

Superscripts 
R  Radiation 

’  Dummy  variable  of  integration 
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